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We perform the a posteriori error analysis of residual type of a transmission 
problem with sign changing coefficients. According to [B] if the contrast is large 
enough, the continuous problem can be transformed into a coercive one. We further 
show that a similar property holds for the discrete problem for any regular meshes, 
' extending the framework from |6]. The reliability and efficiency of the proposed 

estimator is confirmed by some numerical tests. 
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1 Introduction 



Recent years have witnessed a growing interest in the study of diffusion problems with a 
sign changing coefficient. These problems appear in several areas of physics, for example 
in electromagnetism [T2| H~5| H~6| [181 [19]. Thus some mathematical investigations have been 
performed and concern existence results [19] and numerical approximations by the finite 
element methods [191 ffl El E] , with some a priori error analyses. But for such problems the 
regularity of the solution may be poor and/or unknown and consequently an a posteriori 
error analysis would be more appropriate. This analysis is the aim of the present paper. 

For continuous Galerkin finite element methods, there now exists a large amount of 
literature on a posteriori error estimations for (positive definite) problems in mechanics or 
electromagnetism. Usually locally defined a posteriori error estimators are designed. We 
refer the reader to the monographs [21 [31 [T71 [21] for a good overview on this topic. 
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Figure 1: The domain f2 

In contrast to the recent paper [6] we will not use quasi-uniform meshes that are not 
realistic for an a posteriori error analysis. That is why we improve their finite element 
analysis in order to allow only regular meshes in Ciarlet's sense [8]. 

The paper is structured as follows: We recall in Section 2 the "diffusion" problem and 
the technique from [6] that allows to establish its well-posedness for sufficiently large con- 
trast. In Section 3, we prove that the discrete approximation is well-posed by introducing 
an ad-hoc discrete lifting operator. The a posteriori error analysis is performed in Section 
4, where upper and lower bounds are obtained. Finally in Section 5 some numerical tests 
are presented that confirm the reliability and efficiency of our estimator. 

Let us finish this introduction with some notations used in the remainder of the paper: 
On D, the L 2 (D)-norm will be denoted by || • \\d- The usual norm and semi-norm of H S (D) 
(s > 0) are denoted by || ■ \\ St D and | • | Sj d, respectively. In the case D = Q, the index Q will 
be omitted. Finally, the notations a < 6 and a ~ b mean the existence of positive constants 
Ci and C*2, which are independent of the mesh size and of the considered quantities a and 
b such that a < C?\) and C\b < a < C26, respectively. In other words, the constants may 
depend on the aspect ratio of the mesh and the diffusion coefficient (see below). 



2 The boundary value problem 

Let f2 be a bounded open domain of M 2 with boundary T. We suppose that f2 is split 
up into two sub-domains Q + and with a Lipschitz boundary that we suppose to be 
polygonal in such a way that 

see Figure [1] for an example. 

We now assume that the diffusion coefficient a belongs to L°°(f2) and is positive (resp. 
negative) on Q + (resp. f2_). Namely there exists e > such that 

0(2) > eo, f° r a. e. a; G (1) 
a(x) < — e , for a. e. x G (2) 

In this situation we consider the following second order boundary value problem with 
Dirichlet boundary conditions: 

f -div (a Vw) = / in Q, , . 

[ u = on T. 
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The variational formulation of fl3]) involves the bilinear form 

B(u,v) = f aVu ■ Vf 

and the Hilbert space 



n 



#1(0) = { u e H\Q) : u = on T}. 

Due to the lack of coercivity of B on Hq(Q) (see [3 |5j |6j), this problem does not 
fit into a standard framework. In [5J [6J, the proposed approach is to use a bijective and 
continuous linear mapping T from H$ (fl) into itself that allows to come back to the coercive 
framework. Namely these authors assume that B(u, Tv) is coercive in the sense that there 
exists a > such that 

B(u, Tu) > a\\u\\l n Vw G H^(Q). (4) 
Hence given / G L 2 (Q), by the Lax-Milgram theorem the problem 

B(u,Tv)= [ fTv VveH^(tt), (5) 
Jn 

has a unique solution u G Hq(Q). Since T is an isomorphism, the original problem 

B(u,v)= [ fv WveH^Q), (6) 
Jn 

has also a unique solution u G -ffo(O). 

In [6], the mapping T is built by using a trace lifting operator 1Z from H^(J1) into 
H\(VlJ), where £ = <90_ PI (90+ is the interface between 0_ and Q + , 

#i (0 ± ) = {n G F x (0 ± ) : it = on «90± \ £}, 

and 

#i /2 (E) = {n, s : u G #1(0-)} = {u^ : it G #+(0+)} 

is the space of the restrictions to £ of functions in H\(VtJ) (or in H^(Q + )). This last space 
may be equipped with the norms 

i/2,± = inf \u\i,n ± - 

p=U| S 



With the help of such a lifting, a possible mapping T is given by (see 



i> + in + , 

— 1>_ + 272-(f + [s) in 0_, 

where v± denotes the restriction of v to Oj-. With this choice, it is shown in Proposition 
3.1 of [6J that flU) holds if 

\B-(K(vf),K(vf))\ 
Kk = sup — r < 1, (7) 

»efl}(n + ) B + (V,V) 



where B±(u, v) = / aV« • Vv. 

Jn± 

For concrete applications, one can make the following particular choice for TZ, that 
we denote by 1Z P : for any <p G H q (E) we define 1Z p (f) = w as the unique solution 
w G of 

Aw = in u> = on £. 

With this choice, one obtains that Kn < 1 if the contrast 

mhiQ_ \a\ 
maxn + a 

is large enough, we refer to Section 3 of [6] for more details. 

Remark 2.1 Note that in £2, EJ/ the authors consider sub- domains Q + and Q- with a 
pseudo-Lipschitz boundary. However the previous arguments from JBj (shortly summarized 
above) are not valid in this case since the space if 1 (f2 + ) equipped with the norm | • \i t n + is 
not complete. 



3 The discrete approximated problem 

Here we consider the following standard Galerkin approximation of our continuous problem. 
We consider a triangulation T of Q, that is a "partition" of Q made of triangles T (closed 
subsets of Cl) whose edges are denoted by e. We assume that this triangulation is regular, 
i.e., for any element T, the ratio /?-t/pt is bounded by a constant a > independent of T 
and of the mesh size h = max^GT hr, where Ht is the diameter of T and pr the diameter 
of its largest inscribed ball. We further assume that T is conforming with the partition 
of Q, i.e., each triangle is assumed to be either included into Q + or into f2_. With each 
edge e of the triangulation, we denote by h e its length and n e a unit normal vector (whose 
orientation can be arbitrary chosen) and the so-called patch u e = U eC rT, the union of 
triangles having e as edge. We similarly associate with each vertex x, a patch u x = U x& tT. 
For a triangle T, stands for the outer unit normal vector of T. £ (resp. M) represents 
the set of edges (resp. vertices) of the triangulation. In the sequel, we need to distinguish 
between edges (or vertices) included into Q or into T, in other words, we set 

£int = {e E £ : e C Q}, 
£ T = { e e £ : e C T}, 

Mint = {x e Af : x e Q}. 
Problem (JB]) is approximated by the continuous finite element space: 

V h = {v h G H%(n) :v hlT E Ft(T), VT G T} , (8) 
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where t is a fixed positive integer and the space P^(T) consists of polynomials of degree at 
most £. 

The Galerkin approximation of problem (|6]) reads now: Find Uh G Vh, such that 

B(u h ,v h )= [ fv h Vv h eV h . (9) 
Jn 

Since there is no reason that the bilinear form would be coercive on Vh, as in [6] we need 
to use a discrete mapping T h from Vh into itself defined by (see [6]) 



^hVh 



v h+ in Q + , 

-v h - + 2K h (v h+ \Y) in fL 



where 7?.^ is a discrete version of the operator 1Z. Here contrary to [B] and in order to avoid 
the use of quasi-uniform meshes (meaningless in an a posteriori error analysis), we take 

n h = i h n, (io) 

where Ih is a sort of Clement interpolation operator [9] and TZ is any trace lifting operator 
from H^ 2 {T) into (see the previous section). More precisely for ifh G -£T/i(E) = 

W|s ■ v h G V h }, we set 

I h TZ((p h ) = 22 a ^ x ^ 
xeAf- 

where M- = Mint H , \ x is the standard hat function (defined by X x G Vh and satisfying 
X x (y) = $xy) and ftjGl are defined by 



a. 



\u x \ 1 J Wx TZ(ip h ) if x G A/i„t n fL 
<£h(x) if a; G Mnt H E, 



where we recall that u x is the patch associated with x, which is simply the support of X x . 
Note that Ih coincides with the Clement interpolation operator Iq\ for the nodes in f2_ 
and only differs on the nodes on S. Indeed let us recall the definition of Ic\R-((ph) (defined 
in a Scott-Zhang manner [20] for the points belonging to E): 

with 

f \ujx\~ 1 L K{tp h ) if x g Nlm n o_, 

Ac = S 

[ \e x \ 1 J ex TZ((f h )da if x G jV int n E with e x = u x D E. 
The definition of I/j aims at ensuring that 

IhR>{<Ph) = on E. 
Let us now prove that 7^ is uniformly bounded. 
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Theorem 3.1 For all h > and (p n G H n (Yli), one has 

\Rh(<Ph)\i,n- ^ ||^/i||i/2,— 

Proof: For the sake of simplicity we make the proof in the case i — 1, the general case is 
treated in the same manner by using modified Clement interpolation operator. 
Since 1Z is bounded from H^iT,) into one has 

I^Mkn- < ll^lliA- (11) 

Hence it suffices to show that 

|(I-4)^ h )|i, n _<||^||i /2 ,_. (12) 

For that purpose, we distinguish the triangles T that have no nodes in Mint H S to the 
other ones: 

1. If T has no nodes in Mint H E, then I n TZ(iph) coincides with Ic{H(<Ph) on T and 
therefore by a standard property of the Clement interpolation operator, we have 

\{I - I h )H{<p H )\l,T = |(/-/ci)^WIl,T< II^WIIl,^, (13) 

where the patch cut is given by cut = [J T"'- 

2. If T has at least one node in Mint H E, by the triangle inequality we may write 

|(J - 4)ft(<^)|l,T < |(J - / C l)^(^)|l,T + |(/ C 1 - h)n{<p h )\l,T. 

For the first term of this right-hand side we can still use (fT3|) and therefore it remains to 
estimate the second term. For that one, we notice that 



Hence 

\(Ia-h)K(<Ph)\i,T< 
Since lZ(ip n ) = fh on E and due to the definition of Iqi, it follows that for x G T fl E, 



\a x - /3 a 



1 / if h da 



Since all norms are equivalent in finite dimensional spaces, we have for all v n G Pi(e x ), 

\v h (x)\ < \e x \- 1,2 h h \\e x . (14) 
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Moreover, 

^ \tPh\l/2,e; (15) 

where here | ■ \x/2,e x means the standard iJ 1 / 2 (e x )-seminorm. Thus Inequalities f|T4|) with 
Vh = <Ph - \e x \~ l f Cx ¥hda and (JT5J) imply that 

All together we have shown that 

(16) 

Taking the sum of the square of (|T3l) and of ( 1161) . we obtain that 

|(/ - h)K{<p h )\l Q _ < \\K{ip h )\\l Q _ + |^|?/ 2 ,E- 
We conclude thanks to (TTTj) and to the fact that 

\<£h\l/2,Y, ^ ||V 9 /i||l/2-- 

■ 

This Theorem and Proposition 4.2 of [6] allow to conclude that ([9]) has a unique solution 
provided that ([7]) holds, in particular if the contrast is large enough. 
Note that the advantage of our construction of IZh is that we no more need the quasi- 
uniform property of the meshes imposed in [6] . 



-1/2 



(p h da 



4 The a posteriori error analysis 

Error estimators can be constructed in many different ways as, for example, using residual 
type error estimators which measure locally the jump of the discrete flux [21]. A different 
method, based on equilibrated fluxes, consists in solving local Neumann boundary value 
problems (2] or in using Raviart-Thomas interpolant [TJ [TOj [131 E] ■ Here since the coercivity 
constant is not explicitly known, we chose the simplest approach of residual type. 
The residual estimators are denoted by 

= ^kr, Vj = Yl ^It, (17) 
TeT TeT 

where the indicators r/^.T and r\j^ are defined by 

Vr,t = hr\\j T + div(aVu h )\\ T , 

VJ,T = Y ^e /2 HI aV ^' n e]lle ; 

when is an approximation of /, for instance 

h = \T\~ l Jj. 
Note that rf RT is meaningful if q\t €= H^ 1,1 (T), for all TeT. 
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4.1 Upper bound 

Theorem 4.1 Assume that a G satisfies (E])-(0|> and that a\ T G W 1,:l (T), /or a// 

T 6 T. Assume further that |7|) holds. Let u G #q(^) fre the unique solution of Problem 
and let Uh be its Galerkin approximation, i.e. Uh G Vh a solution of (0|). Then one has 

l|V(u- «h) || <VR + VJ + osc(/), (18) 

where 

™<f) = (E^iiz-m 

\rer 

Proof: By the coerciveness assumption (pO), we may write 

II V(« - u h ) || 2 < 5(u - ^, T(ii - Uh )). (19) 

But we notice that the Galerkin relation 

B(u-u h ,v h ) = \/v h G V h 

holds. Hence by taking Vh = Ici^{u — Uh), (JT9~j) may be written 

l|V(u - u ft )|| 2 < B{u - u h , (/ - Ici)T{u - u h )). (20) 

Now we apply standard arguments, see for instance [2T] . Namely applying element- wise 
Green's formula and writing for shortness w = (I — Ici)T(u — Uh), we get 



|V(u-Uh)|| 2 < -J2 I div(aV(u-u h ))w 
TeT 

/ [aV(u — u h ) ■ njw da, 



reminding that w = on T. By Cauchy-Schwarz's inequality we directly obtain 

l|V(u - u h )\\ 2 < ^ ||/ + div (aVuh)|| T |M| T 
TeT 

+ E IH aVn ^ ' "JiWIHIe- 

By standard interpolation error estimates, we get 



|V(u - Uh)|| 2 < ( J2 h T\\f + div (^Uh)||| 
TeT 

+ J2 h e \\la\7u h -n}\\iy /2 \T(u-u h )\ hn . 



Since T is an isomorphism, we conclude that 

\\V(u-u h )\\ < ^h 2 T \\f + dw(aVu h )\\ 2 T 
Ter 

+ h e \\l^n h -n\\\i) 1/2 . 
This leads to the conclusion due to the triangle inequality. ■ 

4.2 Lower bound 

The lower bound is fully standard since by a careful reading of the proof of Proposition 
1.5 of [21] . we see that it does not use the positiveness of the diffusion coefficient a. Hence 
we can state the 



Theorem 4.2 Let the assumptions of Theorems 4-1 be satisfied. Assume furthermore that 



a\T is constant for all T G T ■ Then for each element T G T the following estimate holds 

Vr,t + Vj,t < \u - u h \x >0JT + osc(/, uj t ), 

where 

osc(/,w T ) 2 = h T>\\f - frill* ■ 

T'Clj t 

5 Numerical results 
5.1 The polynomial solution 

In order to illustrate our theoretical predictions, this first numerical test consists in vali- 
dating our computations on a simple case, using an uniform refinement process. Let Q be 
the square (— 1, l) 2 , Q + = (0,1) x (—1,1) and = (—1,0) x (—1,1). We assume that 
a = 1 on Q + and a = fi < on Q_. In such a situation we can take 

K(v+)(x, y) = v + (-x, y) V(x, y) G f2_. 

With this choice we see that 

K n = \fJ-\, 

and therefore for < 1, fl3J) holds and Problem ()6]) has a unique solution. We further 
easily check that the corresponding mapping T is an isomorphism since (T) 2 = T. Similarly 
by exchanging the role of Q + and fi_, (j4]) will also hold if > 1. 
Now we take as exact solution 

u(x, y) = /j,x(x + l)(x - l)(y + l)(y- 1) V(x, y) G Q+, 
u(x, y) = x(x + l)(x - l)(y + l)(y- 1) V(x, y) G 
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/ being fixed accordingly. 

Let us recall that Uh is the finite element solution, and set e^luh) = \\u — Uh\\ and 
^H 1 (uh) = \\u — Uh\\i the L 2 and H l errors. Moreover let us define rj(uh) = T]r + rjj the 
estimator and CV^ (resp. CVm) as the experimental convergence rate of the error ei^i^h) 
(resp. eHi(uh)) with respect to the mesh size defined by DoF^ 1 ^ 2 , where the number of 
degrees of freedom is DoF, computed from one line of the table to the following one. 
Computations are performed with fi = —3 using a global mesh refinement process from 
an initial cartesian grid. First, it can be seen from Tabled] that the convergence rate of 
the H 1 error norm is equal to one, as theoretically expected (see J6J). Furthermore the 
convergence rate of the L 2 error norm is 2, which is a consequence of the Aubin-Nitsche trick 
and regularity results for Problem ([3]). Finally, the reliability of the estimator is ensured 
since the ratio in the last column (the so-called effectivity index), converges towards a 
constant close to 6.5. 



k 


DoF 




CV L 2 


em{uh) 


CVm 


em {uh) 


1 


289 


2.37E-02 




5.33E-01 




6.70 


2 


1089 


5.95E-03 


2.08 


2.67E-01 


1.04 


6.59 


3 


4225 


1.49E-03 


2.04 


1.34E-01 


1.02 


6.53 


4 


16641 


3.73E-04 


2.02 


6.68E-02 


1.01 


6.49 


5 


32761 


1.89E-04 


2.01 


4.75E-02 


1.01 


6.48 


6 


90601 


6.79E-05 


2.01 


2.85E-02 


1.00 


6.47 


7 


251001 


2.45E-05 


2.00 


1.71E-02 


1.00 


6.47 



Table 1: The polynomial solution with fi = —3 (uniform refinement). 



5.2 A singular solution 

Here we analyze an example introduced in [7] and precise some results from [7]. The domain 
Q = (—1, l) 2 is decomposed into two sub-domains Q + = (0, 1) x (0, 1), and — Q \ Q + , 
see Figure[TJ As before we take a = 1 on Q + and a = fi < on fi_. According to Section 
3 of [7], Problem ([6]) has a singularity S at (0,0) if fi < — 3 or if ji e (—1/3,0) given in 
polar coordinates by 

S + (r, 0) = r A ( Cl sin(A#) + c 2 sin(A(^ - 6))) for < 9 < |, 
S_(r, 6) = r A (d x sin(A(# - |) + d 2 sin(A(27r - 6))) for | < 6 < 2tt, 



where A G (0, 1) is given by 



2 ( l-fi 

A = — arccos 



it \2\l + n 
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and the constants ci,C2,dx, d,2 are appropriately defined. 

Now we show using the arguments of Section [2] that for — ~ < fi < and < —3, the 
assumption (|4]) holds. As before we define 

( v + {-x,y) V{x,y) G (-1,0) x (0,1), 
n(v + )(x,y) = I v + (-x,-y) y(x,y) G ( — 1,0) x (—1,0), 
[ v+(x,-y) V(z,y) G (0,1) x (-1,0). 

This extension defines an element of H}_(Q^) such that 

7l(v + ) = v + on E. 

Moreover with this choice we have 

\B4n(v),n(v))\ 

sup — t = 3/i, 

t>#0 

and therefore for 

3 1 /x | < 1, 

we deduce that (jlj) holds. 

To exchange the role of Q + and f2_ we define the following extension from f2_ to for 
w_ G let 

7^(«_)(ar, y) = u_(-a;, y) + v-(x, -y) - v-(-x, -y) V(x, y) G 

We readily check that it defines an element of H+(Q + ) such that 

72-(f_) = u_ on E. 

Moreover with this choice we have (using the estimate (a + b + c) 2 < 3 (a 2 + 6 2 + c 2 ) valid 
for all real numbers a,b,c) 

B + (K(v),K(v)) 
SU P r^l \\ 3 /l^l' 

and therefore for 

3/M < 1, 

we deduce that (jl]) holds. 

For this second test, we take as exact solution the singular function u(x,y) = S(x,y) 
for fi = —5 and \x = —100, non-homogeneous Dirichlet boundary conditions on F are fixed 
accordingly. First, with uniform meshes, we obtain the expected convergence rate of order 
A (resp. 2A) for the H l (resp. L 2 ) error norm, see Tables [2] and [3j There, for sufficiently 
fine meshes, we may notice that the effectivity index varies between 1 and 0.6 for \x = —5 
or between 9 and 6 for /x = —100. From these results we can say that the effectivity index 
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k 


DoF 




CV L 2 


em{uh) 




Tl(Uh) 

P TTl (ill. 1 


1 


289 


1.60E-02 




2.84E-01 




2.57 


2 


1089 


8.66E-03 


0.93 


2.10E-01 


0.45 


1.94 


3 


4225 


4.63E-03 


0.92 


1.55E-01 


0.45 


1.46 


4 


16641 


2.47E-03 


0.92 


1.13E-01 


0.45 


1.09 


5 


32761 


1.80E-03 


0.92 


9.69E-02 


0.46 


0.95 


6 


90601 


1.13E-03 


0.92 


7.68E-02 


0.46 


0.76 


7 


251001 


7.08E-04 


0.92 


6.08E-02 


0.46 


0.61 



Table 2: The singular solution, // = —5, A ~ 0.46 (uniform refinement). 



k 


DoF 




CV L 2 


em{uh) 




em (u h ) 


1 


289 


6.12E03 




1.54E-01 




18.77 


2 


1089 


2.59E-03 


1.29 


9.91E-02 


0.66 


15.04 


3 


4225 


1.08E-03 


1.29 


6.35E-02 


0.66 


12.06 


4 


16641 


4.46E-04 


1.29 


4.04E-02 


0.66 


9.66 


5 


32761 


2.88E-04 


1.29 


3.24E-02 


0.66 


8.65 


6 


90601 


1.49E-04 


1.30 


2.32E-02 


0.66 


7.33 


7 


251001 


7.66E-05 


1.30 


1.66E-02 


0.66 


6.21 



Table 3: The singular solution, fi = —100, A ~ 0.66 (uniform refinement). 

depends on fi, this is confirmed by the numerical results obtained by an adaptive algorithm 
(see below). 

Secondly, an adaptive mesh refinement strategy is used based on the estimator rjT = 
Vr,t + Vj,t, the marking procedure 

f]T > 0.5 md,xi]T> 

and a standard refinement procedure with a limitation on the minimal angle. 

For n = — 5 (resp. fi = —100), Table @] (resp. |5]) displays the same quantitative results 
as before. There we see that the effectivity index is around 3 (resp. 34), which is quite 
satisfactory and comparable with results from [UJ [T4]. As before and in these references 
we notice that it deteriorates as the contrast becomes larger. On these tables we also 
remark a convergence order of 0.76 (resp. 1) in the .f/^-norm and mainly the double in the 
L 2 -norm. This yields better orders of convergence as for uniform meshes as expected, the 
case \i = — 5 giving less accurate results due to the high singular behavior of the solution 
(a similar phenomenon occurs in [TT] for instance). 
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k 


DoF 


e L z{uh) 


CV L 2 


em (uh) 


cv m 


em{u h ) 


1 


81 


2.92E-02 




3.79E-01 




3.39 


5 


432 


3.49E-03 


2.54 


1.40E-01 


1.19 


4.18 


7 


1672 


1.25E-03 


1.52 


8.04E-02 


0.82 


4.07 


10 


5136 


4.26E-04 


1.92 


4.90E-02 


0.88 


3.63 


13 


20588 


1.64E-04 


1.37 


3.14E-02 


0.64 


3.32 


18 


80793 


5.50E-05 


1.60 


1.80E-02 


0.81 


3.23 


24 


272923 


2.39E-05 


1.37 


1.17E-02 


0.71 


2.5 



Table 4: The singular solution, ji = —5, A « 0.46 (local refinement). 



k 


DoF 


eu*(u h ) 


CV L 2 


ejfi (u h ) 


CV H i 


em{ u h) 


1 


81 


1.41E-02 




2.35E-01 




23.59 


4 


363 


1.93E-03 


2.65 


8.77E-02 


1.31 


34.86 


7 


1566 


4.94E-04 


1.86 


4.31E-02 


0.97 


33.10 


11 


5981 


1.23E-04 


2.07 


2.15E-02 


1.04 


33.17 


16 


25452 


2.98E-05 


1.96 


1.05E-02 


0.99 


34.65 


24 


106827 


7.36E-06 


1.95 


5.23E-03 


0.97 


33.89 



Table 5: The singular solution, ji = —100, A w 0.66 (local refinement). 
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